Surface-induced cubic anisotropy in nanomagnets 
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I. INTRODUCTION 
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' We investigate the effect of surface anisotropy in a splierical many-spin magnetic nanoparticle. 

P5 ' By computing minor loops, 2D and 3D energyscape, and by investigating the behavior of the net 

^ , magnetization, we show that in the case of a not too strong surface anisotropy the behavior of the 

many-spin particle may be modeled by that of a macro-spin with an effective energy containing a 
' uniaxial and cubic anisotropy terms. This holds for both the transverse and Neel's surface anisotropy 

models. 

^ . PACS numbers: 75.75.+a, 75.10.HK 

o 

B 

> 
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A magnetic nanoparticle exhibits many interesting and challenging novel properties such as exponentially slow 
Ci . relaxation at low temperature and superparamagnetic behavior above some temperature that depends on the particle's 
^ ■ size and its underlying material. The magnetization of a superparamagnetic particle shuttles in a fast motion between 
the various anisotropy-energy minima. The stability of the magnetization against this thermally-activated reversal has 
become a crucial issue in fundamental research as well as in technological applications. Controlling this behavior, in 
view of room temperature applications, requires a fair understanding of the magnetization dynamics at the nanosecond 
time scale. 

There are two competing approaches to the study of the static and dynamic properties of a nanoparticle. i) One-spin 
particle (OSP): a macroscopic approach that models a nanoparticle as a single magnetic moment, assuming coherent 
rotation of all atomic magnetic moments, and is exemplified by the Stoner-Wohlfarth model for statics and Neel-Brown 
model for dynamics UQ. ii) Many-spin particle (MSP): this microscopic approach involves the atomic magnetic 
moment with continuous degrees of freedom as its building block. It allows taking account of the local environment 
inside the particle, including the microscopic interactions and single-site anisotropy 0. This approach becomes 
necessary when dealing with a very small nanoparticle because the spin non-coUinearities induced by strong boundary 
effects and surface anisotropy invalidate the coherent-rotation assumption. However, investigating the dynamics of an 
MSP is a real challenge. Indeed, within this approach one is faced with complex many-body aspects with the inherent 
difficulties related with analysing the energyscape (location of the minima, maxima, and saddle points of the energy 
I potential). This analysis is unavoidable since it is a crucial step in the calculation of the relaxation time and thereby 
in the study of the magnetization stability against thermally-activated reversal. One may then address the question as 
^ , to whether there exist some cases in which the full-fledged theory that has been developed for the OSP approach [see 
'■^ ' 01 and references therein] can still be used to describe an MSP. However, avoiding somehow the spin non-collinearities 
^ ', induced by surface/interface anisotropy means that some price has to be paid. In Ref. Q it was shown that when 
O 1 the surface anisotropy is much smaller than the exchange coupling, and in the absence of core anisotropy, the surface 
[ anisotropy contribution to the particle's energy is of 4*''-order in the net magnetization components and 2^'^-order in 
^ 1 the surface anisotropy constant. This means that the energy of an MSP with relatively weak surface anisotropy can 
' be modeled by that of an OSP whose effective energy contains an additional cubic-anisotropy potential. 
, It is the purpose of this work to extend the (numerical) result of Ref. |^ to include core anisotropy. This is achieved 
■ by computing minor hysteresis loops, net magnetization components as functions of the field, and energyscape in 2 
and 3 dimensions. We find that because of surface anisotropy the minimum defined by the core uniaxial anisotropy 
splits into 4 minima, reminiscent of cubic anisotropy. We then show that the energyscape can be modeled by an 
effective energy of the net particle's magnetization containing a uniaxial- and cubic-anisotropy terms. We also show 
that this result holds for two different models of surface anisotropy (transverse and Neel). 
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II. MODEL AND COMPUTING METHOD 



We consider a spherical particle of M spins cut from a cube of side N (i.e., — 1 atomic spacings) with simple 
cubic lattice structure. Due to the underlying (discrete) lattice structure, the particle thus obtained is not a sphere 
with smooth boundary because its outer shell presents apices, steps, and facets, resulting in many sites with different 
coordination numbers. 

Our model Hamiltonian is the (classical) anisotropic Dirac-Heisenberg model 

M 

7i = - ^ J,jS, ■ - (.9Mi3)H -^81+ T-ian, (1) 

where is the unit spin vector on site i, H the uniform magnetic field, Af the total number of spins (core and surface), 
and Jij{— J > 0) the nearest-neighbor ferromagnetic exchange coupling. Han is the uniaxial single-site anisotropy 
energy 

Wan =-^if»(s.-e,)2, (2) 

i 

with easy axis and constant Ki > 0. If the spin at site i is in the core, the anisotropy axis is taken along the 
reference z axis and Ki = Kc- For surface spins, this axis is along the radial (i.e., transverse to the cluster surface) 
direction and Ki = Kg. In this case, the model in Q for surface spins is called the Transverse surface anisotropy 
(TSA) model. We use the more general model of transverse direction given by the gradient [the vector perpendicular 
to the isotimic surface — constant defining the shape of the particle, e.g. a sphere or an ellipsoid], because in the 
case of a spherical particle the transverse and radial directions coincide, whereas for another geometry such as an 
ellipsoid they do not. 

A more physically appealing microscopic model of surface anisotropy was introduced by Neel with 

wl'"' = f^EE(«^-^'^)'' (3) 

i j = l 

where Zi is the coordination number of site i and e^j = Vij /rij is the unit vector connecting the site i to its nearest 
neighbors. This model is more realistic since the anisotropy at a given site occurs only when the latter loses some 
of its neighbors, i.e. when it is located on the boundary. The model in Q will be referred to as the Neel surface 
anisotropy (NSA) model 0. 

Qualitatively, the NSA model is not quite different from the TSA model. For example, consider a site i sitting on 
a [100] facet, e.g. in the upper most plane normal to the z axis. It has 4 neighbors on that facet and one below it 
along the z axis. From 0, the corresponding energy reads 



(si • e^f + (sj • -e^f + {si ■ Byf + {si ■ -Qyf + (s^ • -e^f 



= Ks [2sl^ + 2sly + si;\ = 2Kg - Kg 



yi 
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where we have used || ||= 1. This implies that if Kg > the easy direction is along ±8^, i.e., normal to the facet, 
and if Kg < the facet becomes an easy plane. Therefore, upon dropping the irrelevant constant, we rewrite the 
above energy as 

which is the same as the TSA in Eq. ^ for the site considered. More generally, averaging the NSA over a surface 
perpendicular to the direction n leads to [see Eqs. (6, 7) in Ref. 

■HfS^ = -K4\n,\sl + \ny\sl + K\sl) (5) 

thus favoring among (x,y,z) the direction closest to the surface normal. This explains the similarity between the 
results obtained with TSA and NSA. Using the components \na\, a — x,y, z, the atomic surface density reads 



/(n) = ma.x{\n^\,\ny\,\nz\} , 



(6) 
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and thereby we have the particular cases 

Ks/Vi, ^ riy = = 

For comparison, taking account of the atomic surface density, the TSA model is described by 

H^SA = (n . s,)V(n) (7) 

and the whole effect comes from the atomic surface density. Quantitatively, in the NSA model the effect is bigger, 
since for the surface cut perpendicular to the grand diagonal of the cubic lattice the anisotropy completely disappears, 
whereas in Q the surface anisotropy is only reduced by the factor 1 / \/3. 

Because we are dealing with an MSP, the energyscape cannot be represented in terms of the coordinates of all spins. 
Instead, we may represent it in terms of the coordinates of the particle's net magnetization. For this purpose, we fix 
the global or net magnetization, m, of the particle in a desired direction mo (Imoj = 1) by using the energy function 
with a Lagrange multiplier A 1^: 

V s 

= n-MX-{m-mo), m = ' (8) 

To minimize .F, we solve the evolution equations 

= - [si X [si X Fj]] , Fi = -dT/dsi 
A = dT/d\ =-N{m- m-o) , (9) 

starting from s,; = mo = m and A = 0, until a stationary state is reached. In this state m = mo and [s^ x F^] = 0, i.e., 
the torque due to the term J\fX- (m — mo) in J-' compensates for the torque acting to rotate the global magnetization 
towards the minimum-energy directions [see discussion in Ref. 5]]. The orientation of the net magnetization is then 
given either in Cartesian coordinates (ma,, mj,, m^) or in spherical coordinates (On^^n)- 




III. RESULTS AND DISCUSSION 



The anisotropy is taken uniaxial in the core and along the z axis. On the surface, we use the TSA model, except 
in Fig. |H1 where the NSA model is used instead. All anisotropy constants are measured with respect to the exchange 
coupling J, so we define the reduced constants, kc = KqI J-, ks = Kg/ J- In all calculations, unless otherwise specified, 
the magnetic field is applied at an angle of 7r/4 with respect to the core easy axis, i.e., Gz [see discussion at the end 
of section llV| . 

To illustrate the issue stated in the abstract and at the end of section we have chosen kc ~ 0.01, fcj — 0.3, 
which in real units correspond to Kc — 1.3 x 10"^^ J/atom and Ks ~ 4 x 10^^^ J/atom. These values do not 
really correspond to some particular nanoparticle material, since for cobalt nanoparticles, for instance, one would 
have kc — 0.0024. However, they are consistent with the fact that i) kc is typically at least two orders of magnitude 
smaller than exchange, as has been estimated by many authors, ii) the surface anisotropy constant is about one order 
of magnitude larger than the core constant. Here we start with a slightly larger constant, but then it is varied [see 
Fig. □ below]. 

In Fig. n we plot (in solid line) the full hysteresis cycle indicated by narrow- headed arrows, solid for the upper 
and dashed for the lower branch. In triangles is the cycle obtained upon decreasing the field from the positive 
saturation value down to —0.4 and then backwards. It is clear that the appearance of minor loops is an indication 
of the existence of several local minima induced by surface anisotropy. For later discussion, in Fig. [21 we present 
the magnetic structures corresponding to the middle plane of the particle where the indices on top of the structures 
correspond to those shown in Fig.Q] In Fig. Owe see that the bigger is the particle, and thereby the smaller is the 
surface contribution, and the narrower is the minor loop, which implies that the non-coUincarities caused by surface 
anisotropy are indeed responsible for the minor loop. The left graph in Fig. ^ shows how the initial energy minimum 
corresponding to the initial state evolves when the field is varied along the descending branch of the hysteresis loop 
in Fig.^ The graph on the right shows the change in the energyscape when the azimuthal angle is varied at the field 
h = —0.1. Likewise, in Fig. |Slwe plot the evolution in field of the Cartesian coordinates of the net magnetization. 
Upon examining Figs. ^ [S] we see that the orientation of the particle's net magnetization varies according to the 
following scenario: it starts in the minimum defined by the positive saturation state, i.e., (0„ 35°, fn = 0). As the 
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FIG. 1: (Color online) Full and minor hysteresis loops of a spherical particle of A/" = 360 spins with kc — 0.01, ks = 0.3. The 
return field is h = -0.4 {h = H/Hk, Hk = 2Kc). 
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FIG. 2: Magnetic structures corresponding to the middle plane of the particle. The parameters are the same as in Fig. The 
indices on top of the structures correspond to field values shown in Fig. 



field decreases towards negative values, this orientation drifts towards the z axis (Qn ^ 14°,</?n = 0), as confirmed 
by the increasing in Fig. |31 (circles). This continues until the field reaches the value ~ —0.235. Note that the 
net magnetization stays in the x — z plane containing the field direction. Then, the net magnetization jumps to the 
direction (6'„ ~ 38°, <~pn = tt). This occurs through a saddle point in the (y9„ direction as can be seen in Fig. 0| (right). 
As the field further decreases, this direction goes towards the x axis, but still with ipn = tt. When we cycle back from 
the field value h = —0.4 across the value —0.235, the direction of m goes through another saddle point and so does 
not jump back to (0„ ~ lA°,tpn = 0) but continues towards the z axis until the field reaches the value h — 0.15 at 
which it jumps to the state marked by symbols within circles in Fig. |S1 In fact, this is the only state in which the net 
magnetization goes out of the x — z plane. From the latter state, the direction of the net magnetization jumps back 
to the initial state. 
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FIG. 4: (Color online) Left: Energy as a function of the polar angle 9n for (/>„ = and different field values. Right: Energy as 
a function of the polar angle On for h = —0.1 and different values of the azimuthal angle 0n. 



It is clear that this scenario should be dependent on the energyscape, and thereby on the particle's size (see Fig.Ol) 
and surface anisotropy (constant and model). Nonetheless, it shows that surface anisotropy induces more minima 
and saddle points in the energyscape, and thereby more paths along which the direction of the net magnetization can 
evolve as the field is varied. 



IV. EFFECTIVE OSP 



The results discussed above suggest that the MSP may be represented by an effective magnetic moment [as is 
done in the so called macro-spin approximation P, which is just the net magnetic moment of the particle in an 
effective anisotropy potential. Indeed, in Q it was shown that when the surface anisotropy is smaller than exchange, 
in the absence of core anisotropy, the surface anisotropy contribution to the particle's energy is of 4*^ order in the 
net magnetization components and second order in the surface anisotropy constant, and is proportional to a surface 
integral. Now, we may ask the question as to what happens in the presence of core anisotropy. The answer is obtained 
upon computing the 2D and 3-D energyscape for the MSP as explained earlier [see Eq. (O et seq]. The results are 
shown in Fig. |H| They represent the 2D and 3D energyscape of the particle at zero field as a function of the direction 
m of the net magnetization. The locus of m (the unit sphere) has been mapped onto a plane using the azimuthal 
equidistant projection: the spherical coordinates On and become the polar coordinates p = On and 0'^ — ipn of the 
plane. This is the projection used for the UN logo. It preserves the aspect of the "north pole" region, i.e. the region 
around On = 0, which is an easy direction for the core anisotropy. The 2D graph on the right shows that the high 
maxima seen in the 3D graph are indeed at On ~ 7r/2 and that the "beans" on the border are the antipodes of the 
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Coordinates of the net magnetization 
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FIG. 5: (Color online) Cartesian coordinates of the net magnetization m computed along the minor hysteresis loop of Fig. 
The rings mark the unique state at which the magnetization goes out of the x — z plane [see text] . 




FIG. 6; (Color online) Left: Plot of ((6„/7r) cos (6„/7r) sin £) for the same parameters as before. Right; Contour plot 
of the energy in the plane of polar coordinates p — 9„ and 6'„ — ip„ . 

minima around 0„ = and at A, B, C. As can be seen from these plots, the effect of the non-collinearities of the 
magnetization on the energy is to spUt the minimum at On = into four minima at On ^ 28° and = 0, ±7r/2,7r. 
These minima are connected by saddle points at (pn = ±7r/4 and ±3tt/4 and the point at On — becomes a small local 
maximum. The four minima exist over a finite range of applied field, although their positions change continuously 
as a function of the field. They can explain the minor loop of Fig. ^ as follows; In the upper branch of the loop 
(structures 1, 2 and 3 in Fig.|21), the magnetization is in the minimum labeled A. At h ^ —0.235, it jumps to minimum 
C (structure 4). It stays in this minimum (structures 4 — 9) until the field is swept up past h ~ 0.15, at which point 
it jumps to the minimum B (structure 10) and then back to minimum A (structures 11 and 12). The only minimum 
with the magnetization out of the x — z plane, namely B [marked by circles in Fig. is reached when the field is 
swept up, but not when it is swept down. This asymmetry is due to the fact that the field is applied in the direction 
^ri = i'/4, — 0, which is closer to the minimum A than to the other three. 

This suggests that the energy of such an MSP can be modeled by that of a one-spin particle with the net magneti- 
zation m in a potential containing a uniaxial and cubic anisotropy terms, i.e., up to a constant, 

£es. = -Kuniml + Kcuhiml + rriy + m^). (10) 

Indeed, Fig. [7| (first panel, kg = 0.1) shows that for very small kg the energyscape of an MSP is well recovered by 
the effective energy in Eq. H10|l . As kg increases [sec middle panel, kg — 0.3], some deviations start to be seen, and for 
relatively larger values of kg a fit with Eq. H10|) is no longer possible. In fact, in this regime strong deviations from 
collinearity develop, especially near maxima and saddle points, as can be seen on the last panel (fcs = 0.5) in Fig. 13 
In fact, as mentioned earlier, in this case the Lagrange-parameter method in Eqs. JH] |5J) fails because the magnetic 
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FIG. 8: (Color online) Energy as a function of the polar angle 9n for a spherical particle with uniaxial anisotropy in the core 
and Neel anisotropy on the surface. = 0,7r/4 and h = G. The other parameters are as in Fig.0 



state of an MSP can no longer be represented by a net magnetization. Next, fitting the iD energy plotted in Fig. 
with expression (|10|l leads to -ftTuni — 0.0049 and ifcub — 0.0036. The uniaxial-anisotropy term stems mainly from 
the core contribution, though with a small contribution from the surface. Indeed, i^uni = {Nc/M) x i^c, where iVc 
is the number of core spins (= 184 here), leading to Kc — 0.0096, which is close to the constant initially used. This 
simple fitting shows that even for the seemingly strong surface anisotropy (fcs — 0.3) and large surface-to- volume ratio 
(~ 0.49 here), the surface contribution to the effective energy remains smaller than the contribution of the uniaxial 
core anisotropy {Kcuh ^ -ft^uni)- 

The results presented hitherto are all for a spherical particle with uniaxial anisotropy in the core and TSA on the 
surface. Fig. |Slis a plot of the 2D energyscape for a spherical particle with uniaxial anisotropy in the core, as before, 
but now with NSA on the surface. It is clear that the cubic-anisotropy features are also seen in the case of NSA, 
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namely that i) the energy minima are not along the directions (0„ = 0, tt) of the core easy axis, these directions having 
become local maxima, as discussed in the case of TSA, and ii) there is a clear dependence on the azimuthal angle 0„. 

The results above agree and complement those of Ref. and 0, where it was shown that for the TSA and NSA 
there exists a (different) critical value of the surface anisotropy constant that separates i) the OSP Stoner-Wohlfarth 
regime of coherent switching and ii) the MSP regime where the strong spin non-coUinearities invalidate the coherent 
mechanism, and the particle can no longer be modeled by an effective OSP. Obviously, for very small surface anisotropy 
the cubic contribution becomes negligible [see the first panel in Fig. [7| for kg — 0.1] and the Stoner-Wohlfarth 
OSP model provides a good approximation to the many-spin particle. Accordingly, some experimental macroscopic 
estimations of the surface anisotropy constant yield, e.g., for cobalt Ks/J — 0.1 for iron Ks/J — 0.06 
for maghemite particles K^/J ~ 0.04 ^Jl- However, one should not forget that this effective constant depends on the 
particle's size, among other parameters such as the material composition, and for, e.g. a diameter of 2 nm we may 
expect higher anisotropics. 

The present results have been obtained for a field directed at an angle of 7r/4 with respect to the core easy axis. In 
Ref. we showed that the effect of varying this angle from, e.g. to 7r/4, was to decrease the abovementioned critical 
value of the surface anisotropy constant. For simplicity, in this work we only consider this typical value as it is half- 
way between the singular values of and 7r/2, which is also useful for particle assemblies with randomly-distributed 
anisotropy easy axes. 

Finally, we would like to mention that more extensive calculations |l2j have shown that similar features are also 
observed for other crystal structures (fee, bcc, etc.). 



V. CONCLUSION 



We have shown that the behavior of a many-spin particle with core and (weak) surface anisotropy can be modeled by 
that of an effective macro-spin particle whose net magnetization experiences an effective potential energy containing a 
2*^^^- and 4*''-order contributions. However, a question remains open to further investigations: how to distinguish this 
surface-induced 4*''-order contribution from the (usually weak) cubic anisotropy found in many materials in addition 
to the 2**'^-order contribution? In Ref. a fit to the 3D Stoner-Wohlfarth astroid leads to a 4*''-order term an order 
of magnitude smaller than the 2*''-order one. 

The results of the present work could be useful for investigating the thermally activated reversal of a small nanopar- 
ticle. Indeed, the effective OSP with the energy H10|) allows us to include surface effects, though in a phenomenological 
manner, while offering a considerable simplification as compared to the initial many-spin system. However, one should 
note that the effective potential contains a quartic term in the net magnetization components, which renders the anal- 
ysis of the energyscape somewhat more involved but still feasible 0|, as opposed to the case of a many-spin particle. 
In fact, for small enough surface anisotropy for which the effective energy UlUf) holds best, the minima are mainly 
defined by the uniaxial anisotropy and the effect of the cubic-anisotropy contribution is then to modify the loci of the 
saddle points and their number. This obviously renders the calculation of the relaxation rate rather involved |l4j . 

In the general case of arbitrary surface anisotropy intensity, one has to tackle the problem of computing the 
relaxation rates for a many-spin system with all the ensuing complexities. Such difficulties are already encountered, 
in this context, in the much simpler situation of two exchange-coupled spins [T5lll6j|. Obviously, the general situation 
can only be dealt with through numerical approaches. 
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